#include "Matrix.h"


// ##############  Vector  methods  ###############

Vector::Vector(long elements,double val): vector<double>(elements,val)
{

}

Vector::Vector(const Vector& another)
{
    (*this)=Vector(another.size());
	
	for(long i=0;i<size();i++)
	    (*this)[i]=another[i];
		
}

Vector::~Vector()
{



}


Vector& Vector::operator=(const Vector& another) 
{
    resize(another.size());
	
	for(long i=0;i<size();++i)
	{
	    (*this)[i]=another[i];
	}  
	
	return (*this);  

}


bool Vector::operator==(const Vector& another) const
{
    if(size()!=another.size())
	    return false;
		
	long i;
	for(i=0;i<size();i++)
	{
	    if(!nearlyEqual((*this)[i],another[i]))
		    return false;
	}

    return true;
}


bool Vector::operator<(const Vector& another) const
{
    long i;
    for(i=0; i<size() and i<another.size(); i++)
	{
	    if((*this)[i]<another[i])
		    return true;
	}
	
	if(i>=size())
	    return true;
	else
	    return false;
}



bool Vector::nearlyEqualTo(const Vector& another) const
{
    if(size()!=another.size())
	    return false;
		
	for(long j=0;j<size();j++)
	{
	    if(!Vector::nearlyEqual((*this)[j],another[j]))
		    return false;
	}
    
	return true;

}


bool Vector::isDominant(long col) const
{
    for(long j=0;j<size();++j)
	{
	    if(j==col)
		    continue;
			
		if(abs((*this)[col])<abs((*this)[j]))
		    return false;
			
	}
	
	return true;
}


bool Vector::nearlyEqual(double num1,double num2) 
{
    if(abs(num1-num2)<TOLERANCE)
	    return true;
		
	return false;
}

double Vector::largestCoord() const
{
	long till=size();
	double largest = (*this)[0];
	
	for(long i=0;i<till;++i)
	{
		if(nearlyEqual(abs((*this)[i]),abs(largest)))
			continue;
			
		if(abs((*this)[i])>abs(largest))
			largest=(*this)[i];
	}
		
	return largest;

}


void Vector::disp() const
{
    cout <<"\nThe vector :\n";
	
	for(long i=0;i<size();++i)
	{
	    cout <<"\n  " <<(*this)[i];
	}
    
	cout <<"\n\n";
}



// ################   Matrix methods  ####################


Matrix::Matrix(long rows,long coloumns,double val):vector<Vector>(rows)
{
    
	for(iterator iter = begin(); iter!=end(); iter++){
	
	    (*iter)=Vector(coloumns,val);
        }
	

}


Matrix::Matrix(const Matrix& source)
{
    (*this) = Matrix(source.numRows(),source.numCols());
	
	long i,j;
	
	for(i=0;i<numRows();i++)
	{
	    ((*this)[i])=Vector(source[i]);
	}
}
 

Matrix::Matrix(const char * fileName)
{
    ifstream in(fileName);
	
	if(!in)
	{
	    cout <<"\nunable to open file : " <<fileName;
		return;
	}
	
	long rows,cols;
	
	in >>rows >>cols;
	
	(*this)=Matrix(rows,cols);
	
	long i,j;
	for(i=0;i<numRows();i++)
	{
	    for(j=0;j<numCols();j++)
		{
		    in >>(*this)[i][j];
		}
	}
	
	in.close();
}
	


Matrix::Matrix(const Vector& v,Matrix::DIR dir)
{
    
    switch(dir)
	{
	    case RM:
		    
			(*this)=Matrix(1,v.size());
			
			for(long i=0;i<numCols();++i)
			{
			    (*this)[0][i]=v[i];
			}
			
			break;  
			
		case CM:
		
		    (*this)=Matrix(v.size(),1);
			
			for(long i=0;i<numRows();++i)
			{
			    (*this)[i][0]=v[i];
			}
			
            break;
			
		default:
		
		    (*this)=Matrix(1,1,0);
			cout <<"\nUnrecognized second argument in Matrix(const Vector& v,Matrix::DIR dir)\n";
	}
	
}
		


Matrix::~Matrix()
{


}


Matrix& Matrix::operator=(const Matrix& another)
{
    resize(another.numRows());
	
	for(long i=0;i<numRows();++i)
	{
	    (*this)[i]=another[i];
	}

	return (*this);
}


bool Matrix::operator==(const Matrix& another) const
{
    if(numRows()!=another.numRows())
	    return false;
		
	for(long i=0;i<numRows();++i)
	{
	    if((*this)[i]!=another[i])
		    return false;
	}

    return true;
}



bool Matrix::save(const char * fileName) const	
{
    ofstream out(fileName);
	
	if(!out)
	{
	    cout <<"\n unable to open file : " <<fileName;
		return false;
	}
	
	out <<numRows() <<"   "  <<numCols() <<"\n\n";
	
	long i,j;
	for(i=0;i<numRows();i++)
	{
	    for(j=0;j<numCols();j++)
		{
		    out <<(*this)[i][j] <<"   ";
		}
		
		out <<"\n";
	}
	
	out <<"\n";
	
	return true;
}



void Matrix::dispRow(long rownum) const
{
    cout <<"\n ";
	cout.precision(10);
	
	/*
    for(Vector::iterator iter=(*this)[rownum].begin(); iter!=(*this)[rownum].end(); iter++)
	{
	    cout <<*iter <<"  ";
	}  creates problem with const */
	
	long till=numCols();
	for(long i=0; i<till; ++i)
	{
		cout <<(*this)[rownum][i] <<"  ";
	}
	
    cout <<"\n";
}


void Matrix::disp() const
{
    cout <<"\n The matrix :\n";
    for(long i=0;i<numRows();++i)
	{
	    dispRow(i);
	}
}


Matrix Matrix::operator+(const Matrix& b) const
{
    if(numRows()!=b.numRows() or numCols()!=b.numCols())
	{
	    cout <<"\n Matrix dimensions not matching\n";
	    return Matrix(1,1,0);
	}
	
	Matrix c(*this);
	
	long i,j;
	for(i=0;i<c.numRows();i++)
	{
	    for(j=0;j<c.numCols();j++)
		{
		    c[i][j]=c[i][j]+b[i][j];
		}
	}
	
	return c;
}





Matrix Matrix::operator-(const Matrix& b) const
{
    if(numRows()!=b.numRows() or numCols()!=b.numCols())
	{
	    cout <<"\n Matrix dimensions not matching\n";
	    return Matrix(1,1,0);
	}
	
	Matrix c(*this);
	
	long i,j;
	for(i=0;i<c.numRows();i++)
	{
	    for(j=0;j<c.numCols();j++)
		{
		    c[i][j]=c[i][j]-b[i][j];
		}
	}
	
	return c;
}



Matrix Matrix::transepose() const
{
    Matrix transeposed(numCols(),numRows());
	
    long i,j;
	for(i=0;i<numRows();i++)
	{
	    for(j=0;j<numCols();j++)
		{
		    transeposed[j][i]= (*this)[i][j];
		}
	}
	
	return transeposed;
}    



Matrix Matrix::operator*(const Matrix& b) const 
{
    if(b.numRows()!=numCols())
	{
	    cout <<"\n not suitable for multiplication\n";
	    return Matrix(1,1,0);
	}
	
	Matrix c(numRows(),b.numCols());
	
	long i,j,k;
	for(i=0;i<c.numRows();i++)
	{
	    for(j=0;j<c.numCols();j++)
		{
	        for(k=0;k<numCols();k++)
			{
			    c[i][j]=c[i][j] + (*this)[i][k] * b[k][j];
			}
		}
	}
	
	return c;
}


Matrix& Matrix::operator*(double scalar)
{
	
	for(long i=0;i<numRows();++i)
	{
	    for(long j=0;j<numCols();++j)
		{
	        (*this)[i][j]=(*this)[i][j]*scalar;
		}
	}
	
	return (*this);

}


Matrix& Matrix::multRow(long rowNumber,double scalar)
{
    for(long j=0;j<numCols();j++)
	{
	    (*this)[rowNumber][j]= scalar * (*this)[rowNumber][j];
	}
	
	return (*this);
}

Matrix& Matrix::divRow(long rowNumber,double scalar)
{
    for(long j=0;j<numCols();j++)
	{
	    (*this)[rowNumber][j]= (*this)[rowNumber][j] / scalar;
	}
	
	return (*this);
}


Matrix& Matrix::swapRows(long first,long second)
{
   	//why isn't it working ??
    //vector<Vector>::swap((*this)[first],(*this)[second]);
	
	Vector t((*this)[first]);
	
	(*this)[first]=(*this)[second];
	(*this)[second]=t;
	
	return (*this);
}


bool Matrix::isSymmetric() const
{
    long i,j;
	
	for(i=0;i<numRows();i++)
	{
	    for(j=0;j<numCols();j++)
		{
		    if(!Vector::nearlyEqual((*this)[i][j],(*this)[j][i]))
			    return false;
		}
	}
	
	return true;
	
}


bool Matrix::nearlyEqualTo(const Matrix& another) const
{
    if(size()!=another.size())
	    return false;
		
    for(long i=0;i<numRows();i++)
	{
	    if(!(*this)[i].nearlyEqualTo(another[i]))
		    return false;
	}
	
	return true;
	
}


Matrix Matrix::concat(const Matrix& b) const
{
    if(numRows()!=b.numRows())
	{
	    cout <<"\nConcatenation not possible\n";
		return Matrix(1,1,0);
	}
	
    Matrix c(numRows(),numCols()+b.numCols());
	
	long i,j,k;
	for(i=0;i<numRows();i++)
	{
	    for(j=0;j<numCols();j++)
		{
		    c[i][j]=(*this)[i][j];
		}
		
		for(k=0;k<b.numCols();k++,j++)
		{
		    c[i][j]=b[i][k];
		}
	}
	
	return c;
}



Matrix& Matrix::addMult(long first,double multiple,long second)
{
    long count1,count2;
	long j;
	for(j=0;j<numCols();j++)
	{
	    (*this)[first][j]=(*this)[first][j] + multiple * (*this)[second][j];
	}
	
	return (*this);
}
	

bool Matrix::hasGoodPivotIn(long rowNumber) const
{
    if(Vector::nearlyEqual((*this)[rowNumber][rowNumber],0))
	    return false;
		
	return true;
}


bool Matrix ::isColZero(long colNumber) {

      if(numCols() < colNumber) return false;
      
      for(int i=0;i<numRows();++i){
      
           if( (*this)[i][colNumber] != 0 )
	     return false;
	     
      }
      
      
    return true;  
      
}

//     #########   SquareMatrix method definitions


SquareMatrix::SquareMatrix(const char * fileName):Matrix(fileName)
{
    if(numRows()!=numCols())
	{
	    cout <<"\nThe file " <<fileName  <<" does not contain a square matrix\n";
		(*this)=SquareMatrix(1,0);
		exit(0);
	}

	Eigenvalues = new Matrix(numRows(),1);
	eigenvalueCounter=0;
}


bool SquareMatrix::findDominantFor(long row, Matrix& B) 
{
    for(long i=row+1;i<numRows();++i)
	{
	    if((*this)[i].isDominant(row))
		{
		    Vector t=(*this)[i];
			(*this)[i]=(*this)[row];
			(*this)[row]=t;
			
			t=B[i];
			B[i]=B[row];
			B[row]=t;
			return true;
		}
		
	}
	
	return false;	

}


bool SquareMatrix::makeDiagonallyDominant(Matrix& B)
{
    double element;
	
	for(long i=0;i<numRows();++i)
	{
	    if((*this)[i].isDominant(i))
		    continue;
			
		if(!findDominantFor(i,B))
	        return false;
	}
	
	return true;
	
	/*
    for(long i=0;i<numRows();i++)
	{
	    element=(*this)[i][i];    
	    for(long j=0;j<numCols();j++)
		{
		    if(j==i) 
			    continue;
			
		    if(abs(element)<abs((*this)[i][j]))
			    return false;
		}
	}
	*/
	
}


bool SquareMatrix::isDiagonal() const
{
	for(long i=0; i<numRows(); ++i)
	{
		for(long j=0;j<numCols();++j)
		{
			if(i==j) continue;
			
			if(!Vector::nearlyEqual((*this)[i][j],0))
				return false;
		}
	}
	
	return true;
}
